warbleR logo

Streamlined acoustic analysis

Grace Smith Vidaurre and Marcelo Araya-Salas
2015-08-28

Here, we present a case study of how the warbleR functions can be used in a workflow (see diagram below), with some tips on data management in R. For more details on the function arguments, input or output please read the documentation for the function in question (e.g.: help(querxc), ?querxc).

This vignette assumes 1) basic understanding of R, such as loading packages after installment and manipulating objects, 2) the package warbleR has been installed, and 3) the use of R within the RStudio environment.

To start, we will work through how to use warbleR to download recordings from Xeno-Canto. If you have your own recordings (and don’t want to learn about how to obtain sound files from Xeno-Canto), you can skip the first section and go to the Filter recordings by visual inspection section.

Prepare for data collection

Set up a working directory

First, it’s always important to have a working directory for your project. Then, you need to tell R you want to work there. For this example, create a working directory on your desktop. You’ll need to change the username to your own before running the code below.

library(warbleR)

# Create a new directory
dir.create(file.path(getwd(),"warbleR_example"))
setwd(file.path(getwd(),"warbleR_example"))

# Check the location of the directory
getwd()

Obtain metadata and recordings from Xeno-Canto

Next, we can query the Xeno-Canto database for a species or genus of interest. The function querxc has two types of output:

  1. Metadata of recordings accessed by the query word: geographic coordinates, recording quality, recorder, type of signal, etc.

  2. Sound files: Sound files (in .mp3 format) are returned if the argument download is set to be TRUE (by default is FALSE).

You can query Xeno-Canto by genus:

# Query Xeno-Canto for all recordings of the hummingbird genus Phaethornis
Phae <- querxc(qword = "Phaethornis", download = FALSE) 
## Obtaining recording list...
## Processing recording information:
## 695 recordings found!
# Find out what kind of metadata we have
names(Phae) 
##  [1] "Recording_ID"      "Genus"             "Specific_epithet" 
##  [4] "Subspecies"        "English_name"      "Recordist"        
##  [7] "Country"           "Locality"          "Latitude"         
## [10] "Longitude"         "Vocalization_type" "Audio_file"       
## [13] "License"           "URL"               "Quality"
View(Phae)

Or you can query by species

# Query Xeno-Canto for all recordings of the species Phaethornis longirostris
Phae.lon <- querxc(qword = "Phaethornis longirostris", download = FALSE) 
View(Phae.lon)

If you’re interested in the geographic spread of the recording locations, you can use the function xcmaps to visualize locations. xcmaps will create an map image file per species in your current directory, or in the plot window of RStudio if img = TRUE. If img = FALSE maps will be displayed in the graphic device.

# Image type (it) default is jpeg but tiff files are often better resolution
# if you can't open tiff files the "it" argument should be set to "jpeg", or simply use the detault
xcmaps(X = Phae, img = TRUE, it = "tiff") 
xcmaps(X = Phae.lon, img = FALSE)

Filter Xeno-Canto recordings for downstream analysis

Filter recordings by signal type

In most cases, you will need to filter the type of signal you want to analyze. You can filter the recordings prior to downloading them from Xeno-Canto by subsetting the metadata. Then you can input the filtered metadata back into querxc to download only the selected recordings.There are many ways to filter data in R. Shown below is one example that can be modified to fit your own data.

Some of the metadata is not quite consistent across recordings, such as type of signal or recording quality. These are characteristics of the recordings that you will need to explore visually with downstream functions before proceeding with data collection and analysis. However, if you are dealing with large number of recordings, removing lowest quality recordings (D quality level) or selecting specific vocalization types such as song or calls is adviced.

Let say that, for this example, we are interested in investigating the microgeographic variation of long-billed hermit (Phaethornis longirostris) songs. So lets look for a site with the highest number of songs.

# Find out number of available recordings
nrow(Phae.lon) 
## [1] 65
# Find out how many types of signal descriptions exist in the Xeno-Canto metadata
levels(Phae.lon$Vocalization_type)
## [1] "300"         "call"        "calls"       "lekking"     "song"       
## [6] "song at lek"
# How many recordings per signal type?
table(Phae.lon$Vocalization_type)
## 
##         300        call       calls     lekking        song song at lek 
##           1           3           1           1          58           1
# There are many levels to the Vocalization_type variable. 
# Some are biologically relevant signals, but most just 
# reflect variation in data entry.

# Luckily, it's very easy to filter the signals we want 
Phae.lon.song <- droplevels(Phae.lon[grep("song", Phae.lon$Vocalization_type, 
                                ignore.case = TRUE),])

#check resulting data frame
str(Phae.lon.song) 
## 'data.frame':    59 obs. of  15 variables:
##  $ Recording_ID     : Factor w/ 59 levels "10426","107762",..: 58 57 56 46 37 35 33 30 16 14 ...
##  $ Genus            : Factor w/ 1 level "Phaethornis": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Specific_epithet : Factor w/ 1 level "longirostris": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Subspecies       : Factor w/ 3 levels "","baroni","cephalus": 1 1 2 1 1 1 1 1 1 1 ...
##  $ English_name     : Factor w/ 1 level "Long-billed Hermit": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Recordist        : Factor w/ 7 levels "Fernand DEROUSSEN",..: 6 6 5 3 3 3 3 3 3 3 ...
##  $ Country          : Factor w/ 4 levels "Costa Rica","Ecuador",..: 1 4 2 1 1 1 1 1 1 1 ...
##  $ Locality         : Factor w/ 11 levels "Bonampak, Chiapas",..: 4 9 7 5 5 5 5 5 5 5 ...
##  $ Latitude         : Factor w/ 11 levels "10.3333333333",..: 2 10 7 3 3 3 3 3 3 3 ...
##  $ Longitude        : Factor w/ 10 levels "-78.4651","-79.567",..: 7 4 2 8 8 8 8 8 8 8 ...
##  $ Vocalization_type: Factor w/ 2 levels "song","song at lek": 1 1 2 1 1 1 1 1 1 1 ...
##  $ Audio_file       : Factor w/ 59 levels "http://www.xeno-canto.org/10426/download",..: 58 57 56 46 37 35 33 30 16 14 ...
##  $ License          : Factor w/ 5 levels "http://creativecommons.org/licenses/by-nc-nd/2.5/",..: 3 3 5 4 4 4 4 4 4 4 ...
##  $ URL              : Factor w/ 59 levels "http://www.xeno-canto.org/10426",..: 58 57 56 46 37 35 33 30 16 14 ...
##  $ Quality          : Factor w/ 4 levels "A","B","C","D": 1 1 1 1 1 1 1 1 1 1 ...
# Now, how many recordings per locatity?
table(Phae.lon.song$Locality)
#in case you want more than one signal type you could do something like this
Phae.lon[grep("song|call", Phae.lon$Vocalization_type,ignore.case = TRUE),]
# Lets focus on the recordings from La Selva Biological Station
Phae.lon.LS <- Phae.lon.song[grep("La Selva Biological Station, Sarapiqui, Heredia", Phae.lon.song$Locality,
                              ignore.case = FALSE),]

#  and only the high quality one
Phae.lon.LS <- Phae.lon.LS[Phae.lon.LS$Quality == "A",]

# We can check if the location coordinates make sense (all recordings should be from a single place in Costa Rica)
#by setting img=FALSE we can display the map in the graphic device
xcmaps(Phae.lon.LS,img = FALSE)

Once you’re sure you want the recordings, proceed by using querxc to download the files. It is also a good idea to save the metadata as .csv files at this point. This data could be useful later during the analyses and will be definitively needed if you aim for a sceintific publication.

# Download sound files
querxc(X = Phae.lon.LS) 

# Save each data frame object as a .csv file 
write.csv(Phae.lon.LS, "Phae_lon.LS.csv", row.names = FALSE)

Convert Xeno-Canto mp3 recordings to wav format

Xeno-Canto maintains recordings in .mp3, as these are compressed and smaller in size. However, we require .wav format for all downstream analyses. Compression from .wav to .mp3 and back involves information losses, but recordings that have undergone this transformation have been successfully used in research (e.g. Medina-Garcia et al. 2015, see references).

To convert .mp3 to .wav, we can use the warbleR function mp32wav, which relies on a underlying function from tuneR. However, this function does not always work and it remains unclear as to why. This bug should be fixed in future versions of tuneR. If RStudio aborts when running mp32wav, use an mp3 to wav converter online, or download the open source software Audacity (available for Mac, Linux and Windows users). We have made the selected .wav files the available for downloading (see next section).

After .mp3 files have been converted, we need to check that the .wav files are not corrupted and can be read into RStudio (some .wav files cannot be read, despite being in .wav format).

# Neither of these functions requires arguments
# But always check you're in the right directory before executing them
# getwd()
mp32wav() 

#you could use checkwavs to check if files can be read
checkwavs() 

# we will downsample the wav files so the following analyses go a bit faster
# Let's create a list of all the recordings in the directory
wavs<-list.files(pattern="wav$")
lapply(wavs,function(x) writeWave(downsample(readWave(x),samp.rate = 22050),
                                  filename = x))

If you were unable to convert mp3s to .wav format you can download the .wav files with this code:

download.file(url = "https://www.dropbox.com/s/wrh5xuidmvn8rno/wavs.RDS?dl=0",
              destfile = "wavs.RDS",method = "wget")

recs<-readRDS(file = "wavs.RDS")

for(i in 1:length(recs))
  writeWave(recs[[i]],filename = names(recs)[i])

Filter recordings by visual inspection

Note: In case you have your own recordings in .wav format and have skipped previous sections, you must specify the location of your sound files prior to running downstream functions.

The function lspec is a useful tool for:

  • vocal repertoire analyses
  • filtering by visual inspection

This is the first time we can visualize the recordings since downloading from Xeno-Canto, and we can make the most of it.

If your research attempts to assess variation at some social or gepgraphic scales, lspec can provide you with important information about how to steer your analysis. If there is an obvious variation in the structure of vocalizations from different groups (e.g. treatments or geographic regions), you could focus your analysis on a visual classification of vocalizations. You can use lspec to your advantage here, printing spectrograms on paper and classifying signal types by hand.

Whether or not you decide to proceed with visual classification, lspec allows you to visually inspect the quality of the recording (e.g. amount of background noise) or the type, number, and completeness of the vocalizations of interest. You can discard the image files and recordings you no longer want to analyze, as this will become very useful for downstream functions.

# Let's first create a subset for playing with arguments based on the list of wav files we created above
sub <- wavs[c(1,3)]

# flim default, but with ovlp = 10 to speed up process a bit and tiff image files (it = "tiff"), tiff files 
# have better quality and are faster to produce
# if you can't open tiff files the "it" argument should be set to "jpeg", or simply use the detault
lspec(flist = sub ,ovlp = 10, it = "tiff")

# We can zoom in on the frequency axis by changing flim and change the number of seconds per 
#row and number of rows
lspec(flist = sub, flim = c(1.5, 11), sxrow = 6, rows = 15 ,ovlp = 10, it = "tiff")

# Once satisfied with the argument settings we can run all files
lspec(flim = c(1.5, 11),ovlp = 10,  sxrow = 6, rows = 15, it = "tiff")

The image files should look like this:

lspec example, make it 700 pxls width

Notice that the images have the sound file name and page number at the top right corner.

Now we should inspect the spectrograms. Throwing away image files that are poor qualtiy at first glance,(e.g. with lots of background noise) will help us in later steps. The spectro for the file with the recording ID 154143 does not look that good so we will delete these image files. Then we can select only the sound files that have an image file in the folder. This looks a little silly for just one recording, but it could be really useful when working with much larger amounts.

# List the image file in the directory
# Change the pattern to "jpeg" is you used a different image format
imgs<-list.files(pattern = ".tiff") 

# If the maps we created previously are still there, you can remove them from the list easily
imgs <- imgs[grep("Map", imgs,invert = TRUE)]

# Extract the Recording IDs of the files for which image files were kept 
kept <- unique(sapply(imgs, function(x){strsplit(x, split = "-", fixed = TRUE)[[1]][3]}, USE.NAMES = FALSE))

# Now we can get rid of sound files with no image file in the working directory 
snds<-list.files(pattern = ".wav", ignore.case = T) 
file.remove(snds[grep(paste(kept,collapse = "|"), snds, invert = TRUE)])

Data Collection

Automatically detect signals with autodetec

We can move on to data collection with the filtered recordings using the functions autodetec and manualoc. Both these functions work best (faster) on shorter recordings, but there are ways to deal with larger recordings (an hour long or more).

Here are some points that will help us tailor autodetec for our use:

  1. autodetec has 2 types of output:
    • data frame with recording name, selection, start and end times. These are temporal coordinates that will be passed on to downstream functions when you want to measure acoustic parameters. Save this output as an object, or it will not be saved in the environment, but rather sent to the console.
    • the second output is a spectrogram per recording, with red dotted lines marking the start and end of each detected signal, saved in your working directory. In cases in which recordings are long, as in our example, it’s preferable to create long spectrograms ls = TRUE. For shorter recordings (>10 s) short spectrograms my work better (ls = FALSE).
  2. Some important detection parameters to fiddle with:
    • threshold controls detection by relative amplitude (%)
    • bp serves as a frequency bandpass filter
    • msmooth controls combination of window length and overlap to smooth signals that have many peaks and would otherwise be detected as multiple signals
    • mindur & maxdur determines the minimum and maximum duration, respectively, of the signals to be detected

To set these parameters we need to have some idea of the frequency range and duration of the signals we want to detect. Sectrograms produced above can help us to figure this out. It seems that Phaenthornis longirostris songs have frequencies between 2 and 10 kHz and durations between 0.05 and 0.5 s.

If you need to detect all or most of the signals within the recording play around with different argument values to increase detection accuracy. It may be necessary to do several rounds of optimization with different subsets of your recordings. If just a few signals are needed per recordings a low-accuracy detection could still enough selections.

Finally, although autodetec performs automatic signal detection, it doesn’t remove all manual labor from your data collection. Take the time to visually inspect the selections.

# Select a subset of the recordings
wavs <- list.files(pattern = ".wav",ignore.case = TRUE)

# first set a seed so we all have the same results
set.seed(1)
sub <- wavs[sample(1:length(wavs),3)]

# Run autodetec() on subset of recordings, inspect visually until satisfied
autodetec(X=X, bp=c(1,10), threshold=0,  mindur=0.05, maxdur=0.5, envt="abs", flist= NULL,
          msmooth=c(300,90), ls = T, res = 100, flim= c(1,12), wl = 300, set =TRUE,
          sxrow = 6, rows = 15,  redo =FALSE, it = "tiff")

#FIX autodetec TO REDO THE ONES WITH set=T

The image files should look like this (this is the one for the recording ID 154161):

autodetec image example, make it 700 pxls width

Notice that the images have the value of the arguments embedded in their file names. This is useful to try different argument values and compare effect in signal detection.

We won’t save the autodetec ouput in an object until we’re satisfied with the detection. To improve our detection we should play around with the arguments values. In this case we have a good idea of the detection parameters that work best:

autodetec(flist = sub, bp=c(2,9), threshold=20,  mindur=0.09, maxdur=0.22, envt="abs", 
          msmooth=c(900,90), ls = TRUE, res = 100, flim= c(1,12), wl = 300, set =TRUE,
          sxrow = 6, rows = 15,  redo =FALSE, it = "tiff")

This seems to provide a good detection for most recordings (again, this is the one for 154161):

autodetec image example, make it 700 pxls width

Let’s say that we’re satisfied with the detection so we can proceed run the analysis on all the recordings. For this we will remove the argument flist. We should also save the output at this point.

Phae.ad <- autodetec(bp=c(2,9), threshold=20, mindur=0.09, maxdur=0.22, envt="abs", 
          msmooth=c(900,90), ls = TRUE, res = 100, flim= c(1,12), wl = 300, set =TRUE,
          sxrow = 6, rows = 15,  redo =TRUE, it = "tiff")

str(Phae.ad)

#look at the number of selections per sound file 
table(Phae.ad$sound.files)
## 'data.frame':    317 obs. of  5 variables:
##  $ X          : int  1 3 5 7 9 11 13 15 17 19 ...
##  $ sound.files: Factor w/ 6 levels "Phaethornis-longirostris-154070.wav",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ selec      : int  1 3 5 7 9 11 13 15 17 19 ...
##  $ start      : num  0.192 1.466 2.679 3.929 5.624 ...
##  $ end        : num  0.323 1.609 2.838 4.1 5.787 ...
## 
## Phaethornis-longirostris-154070.wav Phaethornis-longirostris-154072.wav 
##                                  17                                 108 
## Phaethornis-longirostris-154123.wav Phaethornis-longirostris-154129.wav 
##                                  30                                  60 
## Phaethornis-longirostris-154138.wav Phaethornis-longirostris-154161.wav 
##                                  22                                  80

Selecting signals based on signal-to-noise ratio (SNR)

Use snrspecs to prepare signal to noise measurements

Filtering your selected signals by signal to noise ratio (SNR) is often a good idea, but not required. Signals that have a ratio close to 1 (or lower) have a very poor quality. A SNR equals to 1 implies that the signal and background noise have the same amplitude.

A SNR filter can be applied at any point in your worklow, after using autodetec or manualoc. However, if you you just need a sample of the signals in each recording, it would make sense to use the SNR functions to perform another quality filter prior to making acoustic measurements. Like the other functions downstream of autodetec or manualoc, the signal to noise functions require the start and end time of the signals.

snrspecs is another function in the family of spectrogram creators. It has very similar arguments to specreator, which we won’t play around with again, but it also has additional arguments for picking a margin over which to measure noise. These margins are very important for calculating SNR, especially when you’re measuring signals with short silence in between. You want to be sure to pick a noise margin that doesn’t overlap neighboring signals.

# A margin that's too large causes other signals to be included in the noise measurement
# Re-initialize X as needed, for either autodetec or manualoc output

# Let's try it on 10% of the selections so it goes a faster
# first set a seed so we all have the same results
set.seed(5)

X <- Phae.ad[sample(1:nrow(Phae.ad),(nrow(Phae.ad)*0.1)),]

snrspecs(X=X, flim = c(2, 11), snrmar = 0.5,mar = 0.7 ,it = "tiff")

The image files should look like this

snrpecs

This margin is far too large! Overlapping the whole signal. Lets try with shorter margins

# This smaller margin is better
snrspecs(X=X, flim = c(2, 11), snrmar = 0.2,mar = 0.7 ,it = "tiff")

snrpecs

It seems that it works. Now we can run the function over all the selections so we can inspect all spectrograms to be sure that the margin(s) you’ve picked are fine

snrspecs(X=Phae.ad, flim = c(2, 11), snrmar = 0.2,mar = 0.7 ,it = "tiff")

Calculate signal to noise ratio for recordings

Once you have picked a margin for all recordings, you can move forward with the SNR calculation. This calculation can allow you to later remove recordings that have a SNR close to 1, as low SNR is indicative of poor quality. Since you’ve already performed several visual filters in the workflow, this step often isn’t necessary, but it can provide you with quantitative information about recording quality.

We will measure SNR on every other selection just to speed up the process

Phae.snr<-sig2noise(X = Phae.ad[seq(1,nrow(Phae.ad),2),], mar = 0.04)

As we just need a few songs to characterize each sound file/individual we could select only the highest SNR selections for each sound file. In this example we will choose the 5 selection with the highest SNR.

Phae.hisnr <- Phae.snr[ave(-Phae.snr$SNR, Phae.snr$sound.files, FUN = rank) <= 5, ]

# Double check the number of selection per sound files 
table(Phae.hisnr$sound.files)
## 
## Phaethornis-longirostris-154070.wav Phaethornis-longirostris-154072.wav 
##                                   5                                   5 
## Phaethornis-longirostris-154123.wav Phaethornis-longirostris-154129.wav 
##                                   5                                   5 
## Phaethornis-longirostris-154138.wav Phaethornis-longirostris-154161.wav 
##                                   5                                   5

At this point would be a good idea to save the selections as a file

#lets save the results
write.csv(Phae.hisnr, "Phae lon autodetec selections.csv")

Manually select signals with manualoc

In some cases manual selection may be preferable, especially if you have shorter recordings or if the automatic detection is not accurate.

manualoc is a function that provides a graphical interface in which you can selec the start and end of the signals. It can often run slowly, depending on the size of the sound files. manualoc can be very useful when you can get away with selecting only a few signals, perhaps if your species has stereotyped signals, or if you’re interested in a specific element of a complex signal across recordings.

We recommend to read the documentation for manualoc prior to running this example. Once you’ve done so, here are some points to keep in mind:

  1. The sole output for this function is a .csv file:
    • contains the time coordinates, selection information and any comments made
    • similar to the autodetec output, these coordinates will be used in downstream functions
  2. Be very precise with your clicks
    • stray clicks will cause manualoc to fail
    • don’t double-click, instead click once and wait for blue bubble to pop up before clicking again
    • suspect a stray click when manualoc starts responding to single clicks
    • if so, stop manualoc with Stop button, or with red Stop button in RStudio console
    • manualoc retains all previous selections in the .csv file and will start up where you left off
  3. Selections can be deleted directly in the manualoc interface
    • use Del-sel button
  4. Selections can also be deleted in the manualoc_output.csv
    • stop manualoc, open the .csv
    • delete the rows corresponding to the unwanted selection(s) prior to starting manualoc again.
  5. Run manualoc within the expected frequency range for your species
    • use argument flim to facilitate signal selection
  6. Run manualoc with oscillograms enabled to improve signal selection
    • when osci = TRUE, the oscillogram or waveform serve as a visual aid
    • use changes in amplitude to select start and end of the signal
    • oscillogram will print to screen when the resolution of the projected spectrogram
      improves (depends on the seltime argument)
    • if seltime = 2, the oscillogram will show up for selections <= 2 seconds

Some other purposes for manualoc:

  1. manualoc can be used in combination with autodetec if you have large recordings:
    • each file will take a long time to load
    • but you can select specific time points to break up the recording
    • then you can feed these time coordinates to autodetec using the data frame argument X
    • this speeds up the automated detection process
    • this can help customize autodetec if you have recordings with different noise or playback treatments
  2. manualoc can also be used for visual classification
  • if you don’t have a printer to print and mark long spectrograms by hand
  • run manualoc with selcomm = TRUE
  • mark individual selections with song or element types using selcomm
  • use specreator to create spectrograms with selcomm text and check visual classifications

Note that you can stop the function at any point by clicking twice on the stop button

# Run manualoc() with frequency range set for Phaethornis longirostris
# recording comments enabled, so you can mark recording quality
# Selection comments enabled to include visual classifications
manualoc(flim = c(1, 11), reccomm = TRUE, selcomm = TRUE, osci = TRUE)

# Read manualoc() output back into RStudio as an object
# This data frame object can be used as input for the functions that follow
manualoc_out <- read.csv("manualoc_output.csv", header = TRUE)

# Note that we don't need to maneuever the data as in autodetec(),
# since we make the desired selections as the function runs. The resulting
# .csv file is equivalent to the output of autodetec() after filtering the 
# automatically detected selections

The graphic device should display something like this

manualoc view

Visualize autodetec or manualoc selections with specreator

specreator serves as yet another option for visual inspection, although at the level of individual selections made through autodetec and manualoc.

Like the other members of the spectrogram-creating family, specreator contains many options related to graphical parameters. With some fiddling around, it’s possible to make images of publication quality. However, some of these graphical parameters do not play well together (especially osci, gr, sc), see the documentation for suggestions.

# Create a subset of 5 recordings analyzed by autodetec() or manualoc()
# Speeds up process of playing around with arguments 
# Run either line below to reinitialize X with either autodetec 
# or manualoc subset as desired

set.seed(50)
X <- Phae.hisnr[sample(1:nrow(Phae.hisnr),5), ]

# Plot selection lines from manualoc() or autodetec()
specreator(X, osci = FALSE, line = TRUE, wl = 300, flim = c(1,11), it = "tiff")

# Change frequency limits of y-axis
specreator(X, flim = c(1, 11), osci = TRUE, line = TRUE, wl = 300, it = "tiff")

# Change width of spectrogram to be proportional to signal duration
specreator(X, flim = c(1, 11), osci = TRUE, line = TRUE, propwidth = TRUE, wl = 300, it = "tiff")

# Change spectrogram size 
# Changing inner.mar and outer.mar arguments improves picsize results
specreator(X, flim = c(1, 11), osci = TRUE, line = TRUE, picsize = 1.5, wl = 300, ovlp = 90,inner.mar = c(4,4.5,2,1), outer.mar = c(4,2,2,1), it = "tiff")

# Run function for all recordings, with final argument settings
specreator(Phae.hisnr, flim = c(1, 11), osci = TRUE, line = TRUE, trel = TRUE, wl = 300,
           ovlp = 90, it = "tiff", res = 300)

specreator example

Measure acoustic parameters

Visualize frequency measurements with trackfreqs

Prior to calculating acoustic measurements, it’s good practice to visualize the accuracy of some important measurements, namely frequency measurements. The function trackfreqs is the last in the family of spectrogram-creators. It allows you to create spectrograms with dominant frequency and fundamental frequency measurements plotted on top of each selected signal.

In general, the fundamental frequency measurements are not as reliable as the dominant frequency measurements. When aocustic measurements are performed in specan, the fundamental frequency reported is the mean of individual fundamental frequency measurements, which is more accurate. Use trackfreqs on all the recordings for which you want to measure acoustic parameters. Scroll through all the spectrograms to get a feeling for how well the frequency measurements will be performed across your recordings.

Like it’s sister functions, trackfreqs has many graphical arguments. It has additional graphical arguments to change colors of the plotting symbols, and size and position of legend labels. These arguments will largely depend on the nature of your selections.

# Note that the dominant frequency measurements are almost always more accurate
trackfreqs(Phae.hisnr, flim = c(1, 11), bp = c(2, 12))

# Play around with the colors and sizes of the symbols
# see par() and points() in RStudio help for more details
trackfreqs(Phae.hisnr, flim = c(1, 11), bp = c(6, 8), col = c("purple", "orange"),
           pch = c(17, 3), res = 300)

# We can change the lower end of bandpass to make the frequency measurements more precise
# If the frequency measurements look acceptable with this bandpass setting,
# that's the setting we should use when running specan() 
trackfreqs(Phae.hisnr, flim = c(1, 11), bp = c(2, 12), col = c("purple", "orange"),
           pch = c(17, 3), res = 300)

trackfreqs image

As fundamental frequency does not seem to be adequately tracked, we will remove acoustic parameters derived from it

Batch-process acoustic measurements with specan

We’re close to finshing the warbleR workflow. We can now perform acoustic measurements with the function specan. This function calculates 22 acoustic parameters across all the specified recordings. It’s a batch process that is much faster than calculating measurements one recording at a time. specan uses and customizes several functions available in the package seewave.

If you require more customized acoustic measurements, you can add your own customizations of seewave functions to the specan code. We suggest that you create your own copy of specan, rather than modifying the package version. Changing the specan code will require getting used to trouble-shooting and debugging in R, unless you’re already comfortable writing your own functions.

# specan() uses the time coordinates in the autodetec or manualoc output
# It will measure acoustic parameters within the start and end times per selection

# Use the bandpass filter to your advantage, to filter out low or high background
# noise before performing measurements
# The amplitude threshold will change the amplitude at which noises are
# detected for measurements 
params <- specan(Phae.hisnr, bp = c(1,11), threshold = 15)

View(params)

str(params)

# As always, it's a good idea to write .csv files to your working directory
## 'data.frame':    30 obs. of  24 variables:
##  $ sound.files: Factor w/ 6 levels "Phaethornis-longirostris-154070.wav",..: 1 1 1 1 1 2 2 2 2 2 ...
##  $ selec      : int  13 15 19 27 29 104 126 132 144 148 ...
##  $ duration   : num  0.176 0.176 0.176 0.159 0.163 ...
##  $ meanfreq   : num  6.53 6.5 6.49 6.61 6.51 ...
##  $ sd         : num  1.77 1.78 1.74 1.84 1.78 ...
##  $ median     : num  6.45 6.34 6.34 6.58 6.42 ...
##  $ Q25        : num  5.51 5.5 5.46 5.39 5.44 ...
##  $ Q75        : num  7.66 7.66 7.61 7.78 7.61 ...
##  $ IQR        : num  2.15 2.16 2.15 2.39 2.17 ...
##  $ skew       : num  1.96 2.26 2.21 2.18 2.04 ...
##  $ kurt       : num  7.73 9.03 9.28 10.05 8.71 ...
##  $ sp.ent     : num  0.943 0.941 0.94 0.947 0.944 ...
##  $ sfm        : num  0.594 0.613 0.572 0.626 0.603 ...
##  $ mode       : num  5.87 6.12 5.88 6.11 6.13 ...
##  $ centroid   : num  6.53 6.5 6.49 6.61 6.51 ...
##  $ peakf      : num  5.84 5.6868 5.7776 6.6079 0.0427 ...
##  $ meanfun    : num  4.04 4.38 4.91 4.81 4.83 ...
##  $ minfun     : num  0.286 0.882 2.756 2.756 2.756 ...
##  $ maxfun     : num  7.35 7.35 7.35 7.35 7.35 7.35 7.35 7.35 7.35 7.35 ...
##  $ meandom    : num  5.9 6.52 6.2 6.87 6.79 ...
##  $ mindom     : num  5.17 5.86 5.34 6.12 5.9 ...
##  $ maxdom     : num  7.19 7.75 7.71 8.01 8.05 ...
##  $ dfrange    : num  2.02 1.89 2.37 1.89 2.15 ...
##  $ modindx    : num  0.447 0.576 0.348 0.627 0.457 ...

Now let’s remove parameters derived from fundamental frequency (based on trackfreqs results)

params<-params[,grep("fun|peakf",colnames(params),invert = TRUE)]

Analysis of geographic variation using specan measurements

Now we can evaluate whether the observed variation in song structure is actually reflected by the acoustic parameters we just measured. For this we will conduct a Principal Component Analysis on scaled (z-transformed) acoustic parameters and look at the grouping of songs (data points) in the scatter plot.

# First run the PCA
pca<-prcomp(x = params[,sapply(params, is.numeric)], scale. = TRUE)

# Check loadings
summary(pca)
## Importance of components:
##                           PC1    PC2     PC3     PC4    PC5    PC6     PC7
## Standard deviation     2.9413 2.2746 1.35134 1.08980 0.9389 0.6195 0.56543
## Proportion of Variance 0.4553 0.2723 0.09611 0.06251 0.0464 0.0202 0.01683
## Cumulative Proportion  0.4553 0.7276 0.82375 0.88626 0.9327 0.9528 0.96968
##                            PC8     PC9    PC10    PC11   PC12    PC13
## Standard deviation     0.52419 0.32400 0.28518 0.24237 0.1744 0.11011
## Proportion of Variance 0.01446 0.00553 0.00428 0.00309 0.0016 0.00064
## Cumulative Proportion  0.98414 0.98967 0.99395 0.99704 0.9986 0.99928
##                           PC14    PC15    PC16      PC17      PC18
## Standard deviation     0.08364 0.06496 0.04952 2.748e-15 3.341e-16
## Proportion of Variance 0.00037 0.00022 0.00013 0.000e+00 0.000e+00
## Cumulative Proportion  0.99965 0.99987 1.00000 1.000e+00 1.000e+00
##                             PC19
## Standard deviation     1.392e-16
## Proportion of Variance 0.000e+00
## Cumulative Proportion  1.000e+00
# Extract PCA scores
pcascor<-as.data.frame(pca[[5]])

# Plot the 2 first PCs
plot(pcascor[,1],pcascor[,2],col=as.numeric(params$sound.files), pch=20, cex=3, xlab = "PC1", ylab = "PC2")

# Add recordings/individuals labels 
x<-tapply(pcascor[,1], params$sound.files, mean)
y<-tapply(pcascor[,2], params$sound.files, mean)

labs <- gsub(".wav","",unique(sapply(as.character(params$sound.files), function(x){strsplit(x, split = "-", fixed = TRUE)[[1]][3]}, USE.NAMES = FALSE)))

text(x, y, labs)

It seems like the songs are grouped by sound file (e.i. individual signature). Now lets look at the song type level. First we need to classified the songs by song type. We can check the spectrograms we previously created to do this.

song types, make it 700 pxl width

Songs from sound files 154070 and 154072 seem to belong to the same song type. Sound files 154129 and 154161 are also from a different song type. Finally, the songs from each of other 2 sound files has a unique structure so each one represent a different song type. We can add this information to the plot by using symbols to represent song types.

# Create a song type variable
# First extract recording ID
songtype<-gsub(".wav","",sapply(as.character(params$sound.files), function(x){strsplit(x, split = "-", 
fixed = TRUE)[[1]][3]}, USE.NAMES = FALSE))

# Now change IDs for letters representing song types
songtype<-gsub("154070|154072", "A", songtype)
songtype<-gsub("154129|154161", "B", songtype)
songtype<-gsub("154123", "C", songtype)
songtype<-gsub("154138", "D", songtype)

# Add song type as a variable representing symbol type
plot(pcascor[,1],pcascor[,2],col=as.numeric(params$sound.files), pch= as.numeric(as.factor(songtype)), 
     cex=3, xlab = "PC1", ylab = "PC2")

# Add song type labels 
x<-tapply(pcascor[,1], songtype, mean)
y<-tapply(pcascor[,2], songtype, mean)

text(x, y, unique(songtype),cex=1.5)

It seems that songs from the same song types are more similar (they are closer together). This also shows that biologically meaningful data can be obtained from sound files that were originally compressed in .mp3 format!

References

Medina‐García, Angela, M. Araya‐Salas, and T. Wright. 2015. Does vocal learning accelerate acoustic diversification? Evolution of contact calls in Neotropical parrots. Journal of Evolutionary Biology. doi: 10.1111/jeb.12694 PDF